Global convection in Earth's mantle : advanced numerical methods and extreme-scale simulations

The thermal convection of rock in Earth's mantle and associated plate tectonics are modeled by nonlinear incompressible Stokes and energy equations. This dissertation focuses on the development of advanced, scalable linear and nonlinear solvers for numerical simulations of realistic instantaneous mantle flow, where we must overcome several computational challenges. The most notable challenges are the severe nonlinearity, heterogeneity, and anisotropy due to the mantle's rheology as well as a wide range of spatial scales and highly localized features. Resolving the crucial small scale features efficiently necessitates adaptive methods, while computational results greatly benefit from a high accuracy per degree of freedom and local mass conservation. Consequently, the discretization of Earth's mantle is carried out by high-order finite elements on aggressively adaptively refined hexahedral meshes with a continuous, nodal velocity approximation and a discontinuous, modal pressure approximation. These velocity--pressure pairings yield optimal asymptotic convergence rates of the finite element approximation to the infinite-dimensional solution with decreasing mesh element size, are inf-sup stable on general, non-conforming hexahedral meshes with "hanging nodes,'' and have the advantage of preserving mass locally at the element level due to the discontinuous pressure. However, because of the difficulties cited above and the desired accuracy, the large implicit systems to be solved are extremely poorly conditioned and sophisticated linear and nonlinear solvers including powerful preconditioning techniques are required. The nonlinear Stokes system is solved using a grid continuation, inexact Newton--Krylov method. We measure the residual of the momentum equation in the H⁻¹-norm for backtracking line search to avoid overly conservative update steps that are significantly reduced from one. The Newton linearization is augmented by a perturbation of a highly nonlinear term in mantle's rheology, resulting in dramatically improved nonlinear convergence. We present a new Schur complement-based Stokes preconditioner, weighted BFBT, that exhibits robust fast convergence for Stokes problems with smooth but highly varying (up to 10 orders of magnitude) viscosities, optimal algorithmic scalability with respect to mesh refinement, and only a mild dependence on the polynomial order of high-order finite element discretizations. In addition, we derive theoretical eigenvalue bounds to prove spectral equivalence of our inverse Schur complement approximation. Finally, we present a parallel hybrid spectral--geometric--algebraic multigrid (HMG) to approximate the inverses of the Stokes system's viscous block and variable-coefficient pressure Poisson operators within weighted BFBT. Building on the parallel scalability of HMG, our Stokes solver demonstrates excellent parallel scalability to 1.6 million CPU cores without sacrificing algorithmic optimality.Computational Science, Engineering, and Mathematic

Advanced Newton Methods for Geodynamical Models of Stokes Flow with Viscoplastic Rheologies

Strain localization and resulting plasticity and failure play an important role in the evolution of the lithosphere. These phenomena are commonly modeled by Stokes flows with viscoplastic rheologies. The nonlinearities of these rheologies make the numerical solution of the resulting systems challenging, and iterative methods often converge slowly or not at all. Yet accurate solutions are critical for representing the physics. Moreover, for some rheology laws, aspects of solvability are still unknown. We study a basic but representative viscoplastic rheology law. The law involves a yield stress that is independent of the dynamic pressure, referred to as von Mises yield criterion. Two commonly used variants, perfect/ideal and composite viscoplasticity, are compared. We derive both variants from energy minimization principles, and we use this perspective to argue when solutions are unique. We propose a new stress-velocity Newton solution algorithm that treats the stress as an independent variable during the Newton linearization but requires solution only of Stokes systems that are of the usual velocity-pressure form. To study different solution algorithms, we implement 2D and 3D finite element discretizations, and we generate Stokes problems with up to 7 orders of magnitude viscosity contrasts, in which compression or tension results in significant nonlinear localization effects. Comparing the performance of the proposed Newton method with the standard Newton method and the Picard fixed-point method, we observe a significant reduction in the number of iterations and improved stability with respect to problem nonlinearity, mesh refinement, and the polynomial order of the discretization.Comment: To appear in Geochemistry, Geophysics, Geosystem

An extreme-scale implicit solver for complex PDEs: highly heterogeneous flow in earth's mantle

Mantle convection is the fundamental physical process within earth's interior responsible for the thermal and geological evolution of the planet, including plate tectonics. The mantle is modeled as a viscous, incompressible, non-Newtonian fluid. The wide range of spatial scales, extreme variability and anisotropy in material properties, and severely nonlinear rheology have made global mantle convection modeling with realistic parameters prohibitive. Here we present a new implicit solver that exhibits optimal algorithmic performance and is capable of extreme scaling for hard PDE problems, such as mantle convection. To maximize accuracy and minimize runtime, the solver incorporates a number of advances, including aggressive multi-octree adaptivity, mixed continuous-discontinuous discretization, arbitrarily-high-order accuracy, hybrid spectral/geometric/algebraic multigrid, and novel Schur-complement preconditioning. These features present enormous challenges for extreme scalability. We demonstrate that---contrary to conventional wisdom---algorithmically optimal implicit solvers can be designed that scale out to 1.5 million cores for severely nonlinear, ill-conditioned, heterogeneous, and anisotropic PDEs

Conformational heterogeneity of the Roc domains in C. tepidum Roc-COR and implications for human LRRK2 Parkinson mutations

Ras of complex proteins (Roc) is a Ras-like GTP binding domain that always occurs in tandem with the C-terminal of Roc (COR) domain, and is found in bacteria, plants and animals. Recently, it has been shown that Roco proteins belong to the family of G-proteins activated by nucleotide-dependent dimerization (GADs). We investigated the RocCOR tandem from the bacteria Chlorobium tepidum with site-directed spin labeling and pulse EPR distance measurements to follow conformational changes during the Roco G-protein cycle. Our results confirm that the COR domains are a stable dimerization device serving as a scaffold for the Roc domains, that in contrast are structurally heterogeneous and dynamic entities. Contrary to other GAD proteins, we observed only minor structural alterations upon binding and hydrolysis of GTP, indicating significant mechanistic variations within this protein class. Mutations in the most prominent member of the Roco family of proteins, leucine-rich repeat kinase 2 (LRRK2), are the most frequent cause of late-onset Parkinson's disease (PD). Using a stable recombinant LRRK2 Roc-COR-Kinase fragment we obtained detailed kinetic data for the G-protein cycle. Our data confirmed that dimerization is essential for efficient GTP hydrolysis, and PD mutations in the Roc domain result in decreased GTPase activity. Previous data have shown that these LRRK2 PD-mutations are located in the interface between Roc and COR. Importantly, analogous mutations in the conserved C. tepidum RocCOR interface significantly influence the structure and nucleotide-induced conformational changes of the Roc domains

Assessing impacts of soil management measures on ecosystem services

Only a few studies have quantified and measured ecosystem services (ES) specifically related to soil. To address this gap, we have developed and applied a methodology to assess changes in ecosystem services, based on measured or estimated soil property changes that were stimulated by soil management measures (e.g., mulching, terracing, no-till). We applied the ES assessment methodology in 16 case study sites across Europe representing a high diversity of soil threats and land use systems. Various prevention and remediation measures were trialled, and the changes in manageable soil and other natural capital properties were measured and quantified. An Excel tool facilitated data collection, calculation of changes in ecosystem services, and visualization of measured short-term changes and estimated long-term changes at plot level and for the wider area. With this methodology, we were able to successfully collect and compare data on the impact of land management on 15 different ecosystem services from 26 different measures. Overall, the results are positive in terms of the impacts of the trialled measures on ecosystem services, with 18 out of 26 measures having no decrease in any service at the plot level. Although methodological challenges remain, the ES assessment was shown to be a comprehensive evaluation of the impacts of the trialled measures, and also served as an input to a stakeholder valuation of ecosystem services at local and sub-national level